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1. Introduction 

The relativistic Fermi gas (RFG) model of Llewellyn-Smith [1] has been used in Monte Carlo 
(MC) generators to describe Charged-Current Quasi-Elastic (CCQE) interactions for the last 40 
years due to its simplicity. The only free model parameter unconstrained by electron-scattering 
data is the axial-mass, Ma, which is well constrained by fits to neutrino-nucleon scattering and pion 
production data as Ma = 1.014 ± 0.014 CeV. Recent measurements of M\ by experiments using a 
heavy nuclear target have found inconsistent results compared with the global dataset [2, 3, 4, 5]. 
In particular, MiniBooNE found a value of Ma = 1.35 ±0.17 GeV with a shape-only fit to their 
data [3], which has led to the name “MiniBooNE large axial-mass anomaly”. 

There has been a great deal of recent theoretical activity on nuclear effects to explain the dis¬ 
crepancy (for a recent review, see reference [6]). There has been a corresponding push within 
the T2K Neutrino Interactions Working Group to develop NEUT [7], T2K’s primary interaction 
generator, to keep up to date with these developments. This report overviews the model develop¬ 
ments in NEUT in Section 2. There is a brief overview of where external data fits enter into the 
general T2K analysis framework in Section 3. A fitting package has been developed to fit NEUT 
MC to published cross-section data. This package is used to fit all of the available CCQE data 
to the new CCQE models available in NEUT as described in Section 4, and the results are used 
to select the default T2K MC model. The Parameter Goodness of Eit (PGoE) test [8] is used to 
test the consistency of the datasets within each model. An error inflation procedure based on the 
PGoE test is described and the final resulfs are presenfed. The besf fil paramefer values and inflaled 
parameter errors will be used as inpufs fo T2K oscillafion fils, and for near detector cross-seclion 
measuremenfs. The work is summarised in Section 5. 

2. NEUT CCQE model developments 

NEUT is Ihe primary inleraclion generator for Ihe T2K experimenl, and developmenfs fo ifs 
inleraclion model musf support all of fhe largel materials used in T2K. The near detector is predom- 
inanfly composed of hydrocarbon and wafer fargels, and Ihe far detector, Super-Kamiokande, is a 
wafer fargef. Additional chemical elemenfs presenf in fhe near defecfor, buf nol used as a primary 
large! for mosf analyses are iron, lead, brass and argon, which musf also be simulaled. 

Two additional nuclear models have been implemented in NEUT to provide allernalives fo fhe 
REG model. The Benhar Speclral Eunclion (SE) model [9] is a more realislic fwo-dimensional 
descriplion of fhe inilial slate nucleon wilhin fhe nucleus in lerms of ifs momenlum and removal 
energy. This SE model includes fhe effecl of short-range correlalions wilhin fhe nucleus, accounting 
for ~20% of fhe lofal cross-seclion. Nole lhal fhe impulse approximalion is slill assumed: fhe 
inleraclion is always wifh a single nucleon. The effective spectral function model of Bodek et 
al. [10] has also been implemented in NEUT, but it is not ready to be a default model for T2K, so 
this is not discussed further in this work. 

Two nuclear effects which are not part of the REG model have also been included. The Ran¬ 
dom Phase Approximation (RPA) is a nuclear screening effect due to long range nucleon-nucleon 
correlations [11]. Two RPA models, relativistic and non-relativistic, are available from the same 
authors. These are different treatments for the quenching of the RPA effect at high four momentum 
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transfer, Q^. The NEUT implementation of both of these RPA models is dependent on and Ey. 
Additionally, a meson exchange current model by Nieves (MEC) [11, 12] has been implemented 
(full details are given in [13]). MEC, or 2p2h models, go beyond the impulse approximation, 
where nucleons within the nucleus are treated as quasi-free. Note that only the lepton kinematics 
are available for the Nieves MEC model, T2K does not have access to the nucleon predictions. An 
effective multi-nucleon ejection model [ 14] is used to simulate the outgoing nucleons, along with 
the NEUT ESI cascade model. As we do not currently have a fully consistent description of the 
hadronic system, we only fit to lepton kinematics in this work. 

With these ingredients, there are two model combinations which can be used as the T2K 
default MC model: 

1. SE-tMEC 

2. REG-tRPA-tMEC 

Note that the SE-i-MEC model is incomplete without RPA. As we do not currently have access to an 
appropriate RPA calculation for the SE model, this incomplete model is the best we have available. 
The BBBA05 form factors [15] are used consistently for both models throughout this work. 


3. T2K analysis structure 

External data fits form an integral part of the T2K analysis structure, as is illustrated by the 
flow diagram in Eigure 1. Eifs fo CCQE and resonanf pion-producfion channels are performed 
independenfly and provide inpufs fo fhe near defecfor (ND280) fif and cross-secfion measuremenfs. 
The ND280 fif consfrains all of fhe flux and cross-secfion pai'amefers which are fhen propagafed fo 
fhe far defecfor for neufrino oscillation measuremenfs [16]. 



Figure 1: A flow diagram illusttating how external data fits are used an an input into various T2K analyses. 


4. CCQE fits to external data 

Eour datasets are used for the CCQE fits presented in these proceedings. These datasets are: 
the MiniBooNE neutrino [3] and antineutrino [17] results on a CH 2 target, which are double¬ 
differential in cosine of the muon angle, cos 6 ^, and the muon kinetic energy, ; and the MINER vA 
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neutrino [18] and antineutrino [19] results on a CH target, which are differential in reconstructed 
four-momentum transfer, 2 q£, defined in Equation 4.1 : 

2(M' -E^ + cos 6^ 

= ->n^, + 2E^^[E^ - ^Ef, - ml cos d^), (4.1) 

where is the muon energy; M„, Mp and nip are the masses of the neutron, proton and muon. 
And where M' = — V, with V as the binding energy of carbon assumed in the analysis. For 

the MINERvA neutrino (antineutrino) dataset, V = 34MeV {V - 30MeV). The restricted phase- 
space MINERvA results, where dp < 20°, were used to reduce model dependence introduced by 
corrected to unsampled regions of phase-space. 

An MC prediction for an arbitrary set of model parameters, x, is produced by reweighting a 
NEUT sample generated using the relevant experimental flux and target. When fitting to the data, 
MINUIT is used to minimise the statistic defined in Equation 4.2: 


N /..DATA /T i' 

k=0 

2 



MiniBooNE v 


MiniBooNE v 


16 16 


-h 


£ £ - vr (5)) (vf - vf^(x)) 

J=0j=0 


MINERvA (4.2) 


where x are the model parameters and is the 16x 16 covariance matrix provided by MINERvA 
which includes cross-correlations between the neutrino (8 bins) and antineutrino (8 bins) datasets. 
Xa and Xp are the normalisation parameters for MiniBooNE neutrino and antineutrino respectively, 
with published normalisation uncertainties of £« (10.7%) and ep (13.0%). Fits can be performed 
to subsets of the data by only including relevant terms in Equation 4.2. When fitting to a single 
MINERvA dataset, the relevant 8x8 portion of V,y is used, cross-correlations between the neutrino 
and antineutrino samples are neglected. Note that the MiniBooNE data has no covariance matrix, 
no bin to bin correlations were published. 

The parameters which can be reweighted in NEUT, and which are investigated as fit parameters 
in this work are: 

• MFC normalisation, as a percentage of the nominal Nieves model. 

• The axial-mass, Ma- 

• The Fermi momentum, pp, which is different for the SF and REG models. 

• The overall normalisation of the CCQE cross-section. 


4.1 Combined fit 

The best fit and parameter values from a combined fit to all four datasets is shown in Table 1 for 
the SF-i-MEC and both RFG-i-RPA-i-MEC models. The relativistic RPA model was favoured over 
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the non-relativistic model, so only the former is considered further in these proceedings. The best 
fit parameter values for the SF+MEC model have an inflated axial mass, are at the limit of no MEC, 
indicating considerable tensions between the parameters favoured by different datasets. None of 
the datasets favour similar parameters when fit with the SE+MEC model individually, which should 
be interpreted as evidence that the SE+MEC model cannot fit all of the datasets in a consistent way, 
despite the comparable obtained for the SE+MEC and relativistic REG+RPA+MEC models. 


Fit type 

X^/DOF 

Ma (GeV) 

MEC (%) 

Pf (MeV) 

7 MB 

Ay 


Rel. RPA 

97.84/228 

1.15+0.03 

27+12 

223+5 

0.79+0.03 

0.78+0.03 

Non-rel. RPA 

117.87/228 

1.07+0.03 

34+12 

225+5 

0.80+0.04 

0.75+0.03 

SF+MEC 

97.46/228 

1.33+0.02 

0 (at limit) 

234+4 

0.81+0.02 

0.86+0.02 


Table 1: Best fit parameter values for combined fits to all four CCQE datasets simultaneously, for the 
RFG+RPA+MEC and SE+MEC fits. 


The best fit distributions for the SE+MEC and relativistic REG+RPA+MEC models are com¬ 
pared with data for MINERvA in Eigure 2 and MiniBooNE in Eigure 3. In the legends of these 
figures, two x^ values are given: the contribution from that dataset to the total, and the total Xmin 
in parentheses. Note that the contribution from each MINERvA dataset ignores correlations, and 
^^MiN total / 2^min V + 2^MiN V’ *^^0 valucs shown in the legends should be treated with some cau¬ 
tion. It is clear from Eigure 3 that MiniBooNE is not dominating the fit, as might be expected given 
the large number of bins in the MiniBooNE datasets. Indeed, the fit exploits the fact that without 
a covariance matrix for MiniBooNE, Xmb ~ 2^min- clear that neither model fits all of the 

datasets perfectly, although this is not reflected by the x^in values in Table 1 . 


xlO^ 



xlO"' 



Figure 2: Comparison between the best fit distributions from the combined fits detailed in Table 1 and the 
MINERvA datasets used in the fit. The x^ value in the legend is the contribution from the dataset shown in 
the histogram, with the total Xmin from all four datasets shown in parentheses. 

4.2 Can these results be trusted? 

A provocative, but important question to ask is whether these fit results can be trusted. It 
has been remarked that the Xmin values suggest very good agreement for both models, but there 
is clearly tension with the data, and the best fit parameter values for the SE+MEC model suggest 
strong tensions between the datasets. 

Without the MiniBooNE covariance matrix, the x^ contributions from MiniBooNE are much 
lower than the number of degrees of freedom it contributes would suggest, and Gaussian statistics 
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(a) Neutrino 



T, (GeV) 


(b) Antineutrino 


Figure 3: Comparison between the best fit distributions from the combined fits detailed in Table 1 and the 
MiniBooNE double-differential datasets used in the fit. The value in the legend is the contribution from 
the dataset shown in the histogram, and the total Xmin from all four datasets shown in parentheses. The solid 
(dashed) lines show the distribution with (without) the normalisation parameters applied. Note that 
some of the cos 0^ slices have been combined in this plot for ease of presentation. 


do not work. This results in two problems: standard goodness of fit tests such as the Pearson 
xlun tsst (SGoF) are unreliable, indicating unrealistically good fits; and hx^ = 1 is no longer 
an appropriate method for calculating parameter errors. We tackle both of these issues with the 
Parameter Goodness of Fit (PGoF) test. 

The PGoF test statistic is defined by Equation 4.3 [8]: 

D D 

ZpGof(x) = zL(x) - £ ^GoF = £ Pr-Pfot, (4.3) 

r=\ 
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where D are the number of datasets, is the minimum x^ in n fit to all D datasets, and Pr and 
Plot are the number of free parameters varied in each fit. 

The PGoF tests the compatibility of different datasets in the framework of the model. Put 
simply, it tests whether fits to subsets of the data pull the best fit parameter far from those found in 
fits to the complete dataset, which would indicate tension between the datasets. A low PGoF value 
implies a poor fit. Although the PGoF test still assumes that the datasets follow a x^ distribution, it 
does not matter how many bins (degrees of freedom) each dataset has. By assuming that without the 
covariance matrix, MiniBooNE data follows a x^ distribution, but with a lower, effective number 
of degrees of freedom, we can use the results of the PGoF test. 

The PGoF results for various subsets of the data are shown for the relativistic RFG+RPA+MEC 
(SF+MEC) model in Table 2 (3). The Xmin value for each subset of the data is found by minimising 
Equation 4.2 with only the relevant terms included. In each fit, p^, MEG normalisation, Ma and 
relevant MiniBooNE normalisation parameters are allowed to vary. An explicit example of how the 
PGoE statistic is calculated is given in Equation 4.4, which corresponds to the final row in Tables 2 
and 3: 

ZpGoF V vs V ~ XmB V + MN V ~ XmB V ~ XmN V (4-4) 



y2. /doe 

SGoF (%) 


PGoF (%) 

All 

97.8/228 

100.00 

17.9/6 

0.66 

MINERvA (V -1- V) 

23.4/13 

3.74 

1.0/3 

79.03 

MiniBooNE (v -t- v) 

58.6/212 

100.00 

2.0/3 

57.69 

V (MB + MIN) 

62.6/142 

100.00 

16.1/3 

0.11 

V (MB + MIN) 

38.5/83 

100.00 

6.1/3 

10.75 

MINERvA vs MiniBooNE 

97.8/228 

100.00 

15.9/3 

0.12 

V vs V 

97.8/228 

100.00 

-3.3/3 

100.00 


Table 2: PGoF statistics for the relativistic RFG+RPA+MEC model. In each fit, M\, MFC normali¬ 
sation and relevant MiniBooNE normalisation parameters are allowed to vary. 

It is clear that there is a great deal of tension in the SF+MEC model. Although there is some 
tension in the relativistic REG+RPA+MEC model, it is clearly a more consistent fit to the data. 
Eor this reason, the SE+MEC model was not selected as the default T2K MC model, and is not 
discussed further in these proceedings. It is reassuring to see good agreement between the neutrino 
and antineutrino datasets for the REG+RPA+MEC model, as shown on the final line of Table 2, 
indicating that the same parameters can be used for neutrino and antineutrino T2K analyses. 

MINUIT uses Ax^ = 1 to define la paramefer errors, which are nof appropriafe for non- 
Gaussian datasets, so the errors returned by our fits are suspicious. This problem has also been 
faced by parton density distribution fitters [20], who overcome it by inflating the value of Ax^ used 
to define fhe eiTors, although there is no general solution offered for choosing that value. 

The PGoE gives a value of the incompatibility between datasets: how much the x^ increases 
between the best fit points of subsets of the data, and the best fit for the combined dataset. The PGoE 
value can therefore be used as a measure of how much the errors must be inflated to cover the differ¬ 
ence between the best fit parameter values from the combined fit, and the best fit values found in fits 
to individual experiments. The rescaled error therefore includes our lack of understanding of the 
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/DOE 

SGoE (%) 


PGoE (%) 

All 

91.SHI'S, 

100.00 

41.1/6 

0.00 

MINERvA (V -1- V) 

12.6/13 

47.75 

1.0/3 

79.49 

MiniBooNE (v -t- v) 

50.2/212 

100.00 

6.5/3 

8.92 

V (MB H- MIN) 

54.8/142 

100.00 

25.1/3 

0.00 

V (MB H- MIN) 

34.1/83 

100.00 

8.5/3 

3.61 

MINERvA vs MiniBooNE 

97.5/228 

100.00 

34.6/3 

0.00 

V vs V 

97.5/228 

100.00 

8.5/3 

3.59 


Table 3: PGoF statistics for the SF+MEC model. In each fit, Ma, Pp^, MEC normalisation and MiniBooNE 
normalisation parameters are allowed to vary. 

differences in the fit results for different subsets of the data. We select Ax^ = Zpgof/^^^pgoF- 
There is some ambiguity over which PGoF value to use, but we use the “MINERvA vs Mini¬ 
BooNE” row of Table 2, as it is the most conservative. This results in a rescaling factor of approx¬ 
imately 2.4. 


4.3 Final fit results 


The correlations between the CCQE parameters for the relativistic REG-i-RPA-i-MEC model 
produced in this work are shown in Table 4. The final errors on the fitted parameters are shown in 
Table 4, where the unsealed errors are also shown for comparison. When an overall normalisation 
parameter for the CCQE interaction channel was included in the REG-i-RPA-i-MEC fit, there was no 
tendency to pull it away from unity. The CCQE normalisation was therefore fixed at the nominal, 
and no error was assigned. 



Figure 4: Correlation matrix for the CCQE parameters for the relativistic REGh-RPAh-MEC model. 


Eit type 

X^/DOE 

Ma (GeV) 

MEC (%) 

Pp (MeV) 

Unsealed 

PGoE scaling 

97.84/228 

1.15±0.03 

1.15±0.06 

27±12 

27±27 

223±5 

223±11 


Table 4: The final, scaled, errors for the CCQE parameters fitted to the relativistic REGh-RPAh-MEC model 
in this work. The unsealed errors are shown for comparison. 

5. Summary 

A number of new models have been implemented in NEUT ready for inclusion in future T2K 
analyses. There are currently two candidate CCQE models which could be used as the default MC 
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model: SF+MEC and RFG+RFG+MEC (with two RPA model options). This work describes fits 
to both models using CCQE data from MINERvA and MiniBooNE to improve the nominal T2K 
model. A generic framework was developed to allow new models and new datasets to be included 
in future iterations of the fit, in order to keep the T2K model up to date. In this iteration, we have 
selected the relativistic RFG+RPA+MEC model as the default CCQE model, and have obtained a 
set of parameters which describe all of the data well, although the best fit parameter values suggest 
that there is some tension between theory and data. 
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